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The multistage factory-warehouse location-allocation 
problem is to decide on locations of warehouses and shipping 
amounts from the factories through the warehouses to meet 
customer demands in such a way that the total fixed plus 
variable costs are minimized. Capacity constraints at facto- 
ries and warehouses are also imposed. 

We present a cost operator algorithm for solving this 
problem. The algorithm takes into account the network struc- 
ture and the submodularity of the objective function. Compu- 
tational results with problems taken from the literature as 
well as new problems are provided. 
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1. Introduction 


In this paper we present a cost operator algorithm 
for solving multistage factory warehouse location— 
allocation problems. The decision variables correspond 
to the warehouse locations and the shipping amounts 
from the factories through the warehouses and into 
the demand centers. The problem is to minimize the 
total fixed costs of locating warehouses plus the total 
variable shipping cost, subject to possible capacity 
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restrictions at the factories and the warehouses, plus 
the demand requirements at the demand centers. 

The algorithm takes advantage of the network struc- 
ture of the supply and demand constraints and the 
submodularity of the objective function. We use cost 
operators [23], to facilitate the movement up and 
down the search tree. This greatly reduces the amount 
of time spent solving transportation subproblems, 
which frequently comprises up to 90% of the compu- 
tational burden. 

In Section 2 we give a problem formulation and we 
discuss some of the more recent research related to 
the multisage location problem. In Section 3 we point 
out lower bounds and fathoming rules obtainable 
from submodular set functions. In Section 4 we trans- 
form the multistage location problem into a transpor- 
tation problem format, and we show how cost opera- 
tors can be used to generate feasible solutions. Sec- 
tion 5 contains a description of the cost operator 
branch and Lound algorithm. In Section 6 we give an 
example, and in Section 7 we provide extensive com- 
putational experience on problems from the litera- 
ture plus some multistage problems of our own. 


2. Problem formulation 


We describe the location problem to be studied in 
this paper as well as similar models which have been 
presented in the literature. We use the following nota- 
tion: 


q  =number of factories; /'= {1, ...,q} 
l = number of warehouses; / = {q + 1, ..., m1}; 
r = number of demand centers; 
m =qtltisna=ltrtls'= {1,8}, 
J= {1+ 1,..., a}. Note that m and 1 have 
each been increased by | so that a dunimy 
factory and a slack demand center can be 
introduced; see Section 4. 
f; = fixed cost of opening warehouse /; 
eee ito ies’ 
warehouse / jes 
cost per unit {factory i to iel' 
Cjj =} of distance (demand centerj /EJ 


shipping from eater ito ie] 
demand centerj /€J 
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A;  =capacity of factory 7; 


S; = capacity of warehouse i: 
d; = demand at location /; 
xj; | = amount shipped from location / to location /; 


1 if warehouse 7 is open, 
0 otherwise. 

The multistage or intermediate location problem, 
which we will call Problem P, can be formulated as: 


Vj = 


Minimize Z(T )= Dy, > ii + 2 Lii ’ (1) 


TCI ieVUT jeJ'v 


subject to 
Xij <A;, i€l'; . (2) 

jes Us 
yy xjj- Sivj SO, (eT: (3) 
jes 

2) xyedj. 1 ES; (4) 
iGI'UT 

Dy Xni~ 21 xy =0 » £67, (S) 
her 
xyj20, i€NUl, jesus; (6) 

Ls. +62, 

y;,= 
“" 10, — otherwise . (7) 


Constrairits (2) and (3) ensure that tiie amount shipped 
out of a factory or a warehouse should not exceed 

its capacity. Constraints (4) require that the demand 
at each demand center be satisfied. Constraints (5) are 
the standard ‘conservation of flow’ constraints which 
require that the amount shipped into each warehouse 
equals the amount shipped out. Constraints (6) and 
(7) are the nonnegativity and integrality constraints 
respectively. The objective in problem P is to minimize 
the total transportation costs from shipping plus the 
total fixed costs of opening warehouses, while satis- 
fying the customer demand. When 


S;< ud; (8) 
jed ; 

we say that warehouse i is capacitated since it cannot 
satisfy all of the demand by itself. The factories can 
also be either capacitated cr uncapacitated depending 
upon the size of A;. Note that in the multistage formu- 
lation when J’ #@, the total warehouse capacity need 
not exceed the total demand since units can also be 
snipped directly from the factories to the demand cen- 
ters. 

The literature on location problems is enormous. 


In this section we discuss recent literature only for 
those problems which are closely related to P. Research 
on multistage location problems has been limited. Geof- 
frion and Graves [12] described and tested an a!gorithm 
for solving multistage multicommodity distribution sys- 
tems using Benders’ decomposition. Ellwein and Gray 
[8] described and tested an algorithm for solving single 
stage (J = Q) location problems and proposed, but did 
not test, an algorithm for solving multistage problems. 
Most of the papers in the literature have concentrated 
on problems which are closely related to the single 
stage version of P, which we denote by P1. Early 
methods attempt to solve PI by imbedding a linear 
programming relaxation of P1 into an implicit enumera- 
tion scheme, see [1,6,7,19,20,21]. Lagrangean relaxa- 
tions have also been used to solve P1, see [4,5,11,16] 
witn some excellent computational results being repor- 
ted by Nauss [17]. 

Recently ‘dual ascent’ methods which utilize pro- 
perties of the dual to various linear programming rela- 
xations of P1 have been used successfully to solve 
uncapacitated location problems [3,9]. Also Guignard 
and Spielberg have generalized the dual ascent methods 
in [3,9] to handle capacitated location problems. 

in the following sections we will describe an algo- 
rithm for solving P which does not solve a relaxation 
of P. Instead we make use of some lower bounds 
obtainable from the submodular property of Z, 
together with some fathoming rules, to enumerate 
implicitly the solutions to P. The movement in the 
search tree from one feasible solution to another is 
facilitated by applying cost operators [22] to P, which 
significantly reduces the amount of time spend in 
solving transportation sub-problems. In the next sec- 
tion we discuss lower bounds on the value of Z which 
are used in the algorithm in Section 5. 


3. Objective function lower bounds 


Many of the lower bound properties which we pre- 
sent here have previously been discussed in (2,10,22], 
in the context of the simple uncapacitated location 
problem. These lower bounds and their properties 
are also useful in solving other kinds of mathematical 
programming problems. 

There are many properties which can be used to 
define submodular set functions. For a discussion of 
them see [10] or [18]. To define one such function 
let Z be a real valued function defined on the finite 
set of subsets of J. 
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For notational convenience we define: 
AZR(A)=Z(A VUB)—Z(A), ACL BCI. (9) 


Then any of the following properties can be used to 
define a submodular function: 


(i) Z(A UB) +Z(ANB)2>Z(A) + Z(B) for all 
subsets 4, BCI. 
(ii) AZ{A) < AZ(B) VA CBCI, i€1~B. (10) 


(iii) Let {B,, .... B,} be a partition of B — A; then 


r 
Z(A) > Z(B) — 27 AZp,(B - By), ACBCI. 
k=] 


(11) 
In most expositions property (i) is taken to be the 
definition of submodularity, and properties (ii) and 
(iii) are shown to be equivalent. Details are omitted. 

In the context of this paper, Z(A ) represents (1). 
The fact that Z(A) is submodular for assignment 
problems was proved by Shapley [21]. (Shapley 
remarks that the extension to transportation prob- 
lems is straightforward.) Note that in P, Z(@) repre- 
sents the value of the solution where all of the ship- 
ments originate from the factories. When /' = @, as in 
Pi, one must be careful in defining 2(A) if A hap- 
pens to be an infeasible set of warehouses, i.e., when 
Lic aS; < Ley dj. In this case we define a "dummy 
factory’ which is always available to service the 
demand centers at a high shipping cost. This makes 
Z(A) large enough so that the value of Z(A ) for any 
infeasible A is at least as large as the value of any 
feasible solution. The dummy factory approach pre- 
serves the submodularity of Z and also yields a dif- 
ferent value of Z for sotitions having different degrees 
of infeasibility. As we will see later, this is helpful in 
deciding which warehouses to open when we are 
working with infeasible sets of warehouses. 

Properties (10) and (i 1) can be interpreted as 
adding or subtracting warehouse from a given set A 
of open warehouses. For example property (10) says 
that the addition of warehouse i to the set A decreases 
the total cost by at least as much as the addition of 
warehouse i to the set B when A C B. Thus it is similar 
to the ‘decreasing returns to scale’ condition in eco- 
nomics. 

Properties (10) and (11) can be used to characterize 
solutions to P and to derive lower bounds on the value 
of any solution to P. The algorithm to be described in 
Section 5 searches for sets T C / which have the follo- 
wing two properties: 


(i) AZ(T)2>0, ViEl-T, 
(ii) AZ(T - fi})<0, VIET. 


Seis T C J which satisfy (i) and (ii) are such that the 
addition to 7 or the deletion from T of a single facility 
does not cause Z to decrease. Clearly any optimal solu- 
tion to P satisfies (i) and (ii). In a sense properties (i) 
and (ii) characterize the set of all ‘locally optimal’ solu- 
tions to P. The globally optimal solution is the vest 
locai solution, which must be found by a search pro- 
cess. One of the factors which make it hard to find 
the globally optimal solution is that there are many 
locally optimal solutions which have nearly optimal 
objective function values. Thus in any enumeration 
procedure, a considerable amount of effort is normally 
required to eliminate these nearly optimal solutions 
from consideration. In a practical sense however, it 
may be true that the nearly optimal solutions, say 
those within 1% of optimality, may indeed be as 
valuable or ‘as optimal’ to a decision maker as a glo- 
bally optimal solution. Given the inaccuracies in the 
cost data and other environmental and political factors 
which must be taken into consideration, it would be 
desirable to have not only an optimal solution to P, 
but in addition a list of solutions which are nearly opti- 
mal. We will point out how this may be accomplished 
using the algorithm of Section 5. 

In the initialization phase of any location algorithm, 
two rules can be applied in order to permanently open 
or close warehouses. 


Rule 1. /f AZ, — {i}) <0 forall i EI” for some sub- 
set I* Cl, then there exists an optimal solution with 
all warehouses i € I” open. 


Proof. Let T be an optimal set of warehouses and 
assume that T-/* = {i;,..., ig} #@. Define Ty = T, 
and T, =7;_, U {ix}, & = 1, ...,q. Then by property 
(10) we have 


AZ (Tx 1) = ZT) — 2(Te 1) $ AZ; ~ {ie SO, 
k=1,..5q 

and thus TU /" is also an optimal set of warehouses. 

Rule 2. [f AZ,(/*) > 0 for alli€C* forsome C* Cc 


1 — 1", there exists an optimal solution with all ware 
houses i € C* closed. 


The justification for Rule 2 is similar to that for 
Rule 1. We will show, in the next section how the tes- 
ting of Rule 1 and Rule 2 requires only the applica- 
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tion of two cost operators for each warehouse. In 
many cases, as will be seen in Section 7, Rules 1 and 2 
can be used to fix open and closed a large portion of 
the warehouses in an optimal solution. 

Another property of submodular set functions 
which can be derived from (10) and (11) is the follo- 
wing: Let {A,,...,.A,} and {Qj, ..., Q,} be partitions 
of A — T such that Q; C A; for each i= 1, ..., ¢ and 
some j= 1, ..., 7. Then 


UT)> ZA) ~ 24 AZ4,(4 ~ Ax) 


t 
> Z(A) — 2 AZo, (A ~ Ox) (12) 


for all subsets 7 satisfying TC A C /. The quantity 
on the right-hand side of (11) provides a value of 
Z(T), for any T C J. Property (12) says when we use 
a more refined partition of A — T, we get a weaker 
lower bound. Suppose we let A = /. Then the most 
refined partition of J — Tis clearly {{i,}, ..., {i,t}, 
where i, © / — Tand t = |/ — T|. Then for this parti- 
tion (12) yields: 


AT) > ZN) - 2 AZ: ~ tigh). (13) 


Among the class of lower bounds (11), (13) is the 
weakest. Nutice that once Rule | has been applied 

all of une terms on the right-hand side of (13) have 
oeen calculated. Thus (13) can be applied at any time 
after Rule 1 has been tested without any extra compu- 
tational effort. A stronger lower bound than (13) 
could be obtained with some added computational 
effort by partitioning / into sets, A;,, of size two and 
calculating AZ 4, (7 — Ax). We have not yet tested this 
idea. 

Another fathoming device which provides an upper 
bound on the maximum number of warehouses in an 
optimal solution to P can be obtained as follows: Let 
r=|J—/" |). Let fx, Sus olky k;€1—1" , be a nonde- 
creasing CLESnE}O! of the fixed chintges not in/”. 

Define 


DV;, = AZ; {I oF {i}) — Si, ’ 
and let DV), > + ,» DV;, be nonincreasing ordering of the 
DV;. “(Note that DV;, <0.) The quantity —DV,,, repre- 
sents the smallest sosible incremental savings on the 
total shipping cost when warehouse /; is opened. Let 
Z" be any upper bound on the optimal objective func- 


Rel-0" (14) 


tion value for P. Then an upper bound on the number 
of warehouses in an optimal solution to P (contain 
I*)is [1"| +k", where 


k* t—k 
2) — 24 f+ Qu fet 22 IDV,1< 
1el [=1 sz 


"+1 t—k*-1 


<2" <2) Lt Le fig * a IDV;,I - 
(15) 


Note that the first two terms on the right-hand side 
of (15) give the smallest possible shipping cost, the 
next term is the smallest sum of k* + 1 fixed charges, 
and the last term is the smallest possible cost of closing 
down t -- (k* + 1) warehouses. Therefore the right- 
hand side of (15) is a lower bound on the value of an 
optimal solution to P containing |/* | +k* + 1 ware- 
houses. If this lower bound exceeds Z", then ariy 
optimal solution will contain at most |/*|+k* open 
warehouses. 

In the algorithm of Section 5 we use as fathoming 
devices bounds obtained from (3) on the value of 
an optimal solution, and bounds obtained from (15) 
on the number of warehouses in an optimal solution. 


4. Use of cost operators to solve problem P 


Given the choice of a subset TC J, problem P 
becomes an ordinary transshipment problem. Rather 
than resolving this problem each time T changes, we 
use the operator theory of parametric programming 
[23] to change the problem and derive the new opti- 
rnal solution simultaneously. 

We shall say that cell (i, /) has been fixed out of the 
basis when a cell cost operator has been applied to the 
problem and its solution so that the cost c,; has been 
driven to +M, where M is so large that x,; = 0 in any 
optimal solution to the new problem. 

We also say that cell (i, /) has been fixed in the basis 
when a cell cost operator has been applied to the prob- 
lem and its solution so that the cost ¢;; has been driven 
to —M, where M is so large that x;; = Min(S;, d;) in any 
optimal solution to the new problem. 

Fig. 1 shows an example of a 2 X3 X4 multistage 
location problem. Rows 1 and 2 correspond to the 
factories, rows 3, 4, and 5 to the warehouses, and 
row 6 represents the dummy factory. Columns 1, 2 
and 3 correspond to the warehouses, columns 4 through 
7 to the demand centers, and column 8 is a slack 
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Warehouses Demand Centers 


Factories 


Ware- 
houses 


Dummy 
Factory 


Fig. 1. 


column. Rows 1 and 2 contain the costs of shipping 
from the factories to the warehouses and from the fac- 
tories directly to the demand centers. Cells (3, 1), 

(4, 2), and (5, 3) contain the fixed charge of the cor- 
responding warehouse, divided by its capacity. In any 
solution to P, these cells will contain the unused ware- 
house capacity. Cells (3, 2), (3, 3), (4, 1), (4, 3), G, 1) 
and (5, 2) cannot be used because of their large costs; 
in effect these arcs have been removed from the prob- 
lem. The cells in rows 3 through 5, and columns 4 
through 7, contain the costs of shipping from each of 
the warehouses to each of the demand centers plus 
the proportional fixed charges. Notice that the cells 

in dummy factory row 6 and slack column 8 contain 
costs equal to M or —M. This can be explained in the 
following manner. To obtain the solution, where: 


y,=0 set Ci, = —M and Cm i-g = —M; 
y¥,= 1 Set Cin = M and Cmj-g =M 


(in Fig. 1, m= 6,n =8,q=2,andi€ {3, 4, 5}). 

By solving the transportation problem in Fig. 1 we 
would obtain an optimal solution to P when T= {2, 3}, 
i.e., when warehouse 1 is closed and warehouses 2 and 
3 are open. Notice that cell (3, 8) has a cost of —M, 
thus in the optimal solution x33 = 53, which has the 
effect of closing down warehouse. Cell (6, 1) has a 
cost of —M which causes all of the demand in column 
1 to be satisfied by the dummy factory, i.e., x6, = 53. 
Thus setting c33 = —M and cg, = —M effectively closes 
down warehouse 1. On the other hand for warehouse 


2,¢4g = Mand c62 = M, so that xqg = X62 = 0. This 
causes the demand in column 2 to be satisfied from 
rows 1, 2, or 4. The shipments from rows | and 2 
represent units being shipped from factories | and 2 

to warehouse 2. The shipment from row 4 in column 

2 represents the unused warehouse capacity, and it 

is charged at the proportional tixed charge rate. Notice 
that if cell (4, 2) contains x42 units then exactly Sq — 
X42 units can be used to ship from warehouse 2 to the 
demand centers. This is exactly the amount which is 
made available to warehouse 2 from the factories. Thus 
the conservation of flow equations (5) have been satis- 
fied. Because a proportional amount of the fixed charge 
is assigned to both the used and unused parts of the 
warehouse capacity, the total fixed charge is covered 

in any feasible solution. 

Fig. 4 contains an optimal solution to the 2 X5 X4 
examples solved in Section 6. In Fig. 4, T= {1,3} so 
that warehouses | and 3 are opened, and warehouses 
2,4 and 5 are closed. Factory | is not used at all so 
that x; .:9 = 61. Factory 2 ships 21 units to warehouse 
1, 40 units to warehouse 3, and 16 units directly to 
demand center 1. War-house 1 ships all 21 of the 
units it received from Factory 2 to demand center 3. 
The remaining unused 4 units of capacity at ware- 
house | are in cell (3, 1) and are charged at the propor- 
tional fixed charge rate. All of the flow for warehouses 
2,4, and 5 isin the last column, and their demand is 
satisfied entirely by the dummy factory. Warehouse 3 
is used to capacity shipping 22 units to demand center 
2, and 18 units to demand center 4. In this example 
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Warehouses 


Factories 


Ware- 
houses 


Dummy 
Factory 


Fig. 2. 


each of the demand centers is supplied from a single 
warehouse. This is not the case in general. 

The idea behind the cost operator approach to 
solving P(T) is to open the warehouses in T and close 
those in / — T by fixing in or out of the basis each of 
the cells (m, 1) through (m, J) and (q + 1, 2) through 
q+l,n). For example the problem P(/), in which all 
of the warehouses are open, can be obtained by sol- 
ving the transportation problem shown in Fig. 2. 
Then to apply Rule 1 in Section 3 we use two cost 
operators for each i € J. For example we would calcu- 
late AZ,(7 — {3}) = Z() — ZU — {3}) by fixing in 
cells (3, 8) and (6, 1). Doing this would yield a solu- 
tion to the problem shown in Fig. 1, and permits the 
evaluation of Rule 1. In general, the addition of a 
warehouse to a given set or the deletion of a ware- 
house from a given set requires the fixing in or out 
of two cells. The amount of work needed to fix two 
cells in or out using cost operators is much less than 
the computational effort of solving a transshipment 
froblem from scratch. Since the branch and bound 
search algorithm to be described in the next section 
requires the solution of transshipment problems for 
many sets 7, the total computational effort saved by 
the fixing in and fixing out procedure is very large. 

Also we should mention that in the case where 
1’ =, as in single stage problems, only one cost opera- 


tor is needed to open or close a warehouse. In the single 


stage formulation columns 1 through / and rows 1 
through 2 are absent. Thus in order to fix in or out 
warehouse / we need only apply a cost operator to 
cell (i, n). 


Demand Centers 


5. The cost operator algorithm 


To describe the cost operator branch and bound 
algorithm we use the following notation: 


l = level of the search tree; 

Q = =set of open warehouses; 

List(/)= list of warehouses which may be opened 
on level J; . 

Z“  =current upper bound; 

X ~~ =current best solution. 


We begin by setting / = 1,Q=@ and Z" = Z(@) 
which is the (high) cost of supplying all demar:ds from 
the factories directly. Then Rule 1 (see Section 3) is 
applied and we let Q=/". Then for each iG / — Q we 
apply the lower bound test (13) by setting T= QU 
{i}, and checking to see whether 


Z° >Z(I)— > AZ,(I — {k}). (16) 


If the right-hand side of (16) (which is 2 lower bound 
on the value of Z(Q U {i})) exceeds Z", then ware- 
house i can be permanently closed. For each i for 
which (16) holds true we apply cost operators to cal- 
culate AZ,(Q). If AZ,(Q) > 0 then we permanently 
close warehouse i. Then we let List(1) be the set of 

all i such that AZ,(Q) <0. (We place the warehouses 
in List(1) in order of non-increasing objective function 
values.) Next we remove the last warehouse, i, from 
List(1) and replace Q by Q U {i}. We let Z" = Z(Q), 
update X. If List(1) is empty, then we stop. Other- 
wise we calculate k* as in (15). If k* = 1, then we stop. 
Otherwise we replace / by / + 1 and continue to the 
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next level. At each level after the first, we perform 
the lower bound test (16) for each warehouse, i, in 
List(/ — 1) by letting T= Q U {i}. If warehouse / passes 
the test (16), we calculate AZ,(Q). Again we let 
List(Z) be the set of warehouses such that AZ,(Q) <0 
in nonincreasing order. If List(/) is not empty, we 
remove i, the last element of List(/); let ZY = Z(Q U 
{i}) and update X. If List(/) is now empty, we back- 
track by replacing / by / — 1, removing the last ware- 
house placed in Q and then continuing as before. 
Otherwise we calculate k* as in (15). If />k*, we 
backtrack as just described. Otherwise we replace 
lby / + 1; replace Q by Q U {i}, and continue as 
before. The algorithm terminates when List(1) has 
been emptied. 

Now we formally state the cost operator multistage 
location—allocation algorithm 


Step 1. Initialization. Set |= 1;z" = Z(Q). Apply Rule 
1. Let Q =/". For each i €/ — Q apply test (16) let- 
ting T=Q U {i}. For each i which satisfies (16) cal- 
culate AZ,(Q). Order those i such that AZ,(Q) <0 in 
nonincreasing order and place them in List(1). If 
List(1) is empty, go to Step 5. Otherwise let k be the 
last warehouse in List(1). Let Z" =Z(Q VU {k}); 

update XY. Remove warehouse k from List(1). If List(1) 
is empty, goto Step 5. Otherwise go to Step 2. 


Step 2. Bound by k*. Calculate k* as in (15). If />k", 
go to Step 4. Otherwise replace / by / + 1; replace Q by 
QU {k} and go to Step 3. 


Step 3. Forward branching. For each warehouse i in 
List(/ — 1) apply test (16) letting T= QU {i}. For 

each i which satisfies (16) calculate AZ,(Q). List those 

i such that AZ(Q) <0 in non-increasing order and place 
them in List(/). If List(/) is empty, go to Step 4. Other- 
wise let kK be the last warehouse in List(Z). Let Z° = 

Z(Q U {k}); update X. Remove warehouse k from 
List(/). If List(Q) is empty, go to Step 4. Otherwise go 

to Step 2. 


Step 4. Backtracking. Replace 1 by 1 ~ 1. Remove 
from Q the last warehouse placed in Q. If 7 = 1 and 
List(1) contains only one warehouse go to Step 5. 
Otherwise let k be the last warehouse in List(/). 
Remove k from List(/). If List(J) is empty, go to Step 
4. Otherwise replace Q by Q U {k}; replace? by / + 1 
and go to Step 3. 


Step 5. Termination. Stop. The current solution X, is 
optimal. 


Notice that each time Step 3 is performed a new 
list of feasible solutions to P is available. Many times 
these solutions have nearly optimal objective function 
values and thus they might be worthwhile to save. 

After calculating /*, the cost operator algorithm 
proceeds to open, one at a time, those warehouses 
which cause the largest decrease in the objective func- 
tion value. This is continued until at level /, no further 
decrease in the objective function value is possible (i.e. 
List(/) is empty). Then we move to the backtracking 
Step 4. The best feasible solution obtained by the 
cost operator algorithm before the jirst backtracking 
is called the greedy solution. We denote by Z“ the 
objective function value of the greedy solution. In 
[18], a worst case bound on Z© was derived for sub- 
modular set functions which is, 


7-27 (ie oye 
ZQ)—-Z" \ "+k" 
where Z” is the optimal objective function value and 
k* is obtained from (15). 

In practice however, the actual percentage error 
obtained by the greedy solution is much smaller 
than the worst case bound. In Section 7 we will see 
that in most cases the greedy solution value is within 
0.5 percent of optimality. 

Finally we should point out that it would be trivial 
to ‘reverse’ the cost operator algorithm so that 
instead of starting with all of the warehouses closed, 
and opening them one at a time, we could start with 
all of the warehouses open and close them one at a 
time. This could be done by defining Z'(T) = Z(/ -- T) 
in (1). Z'(T) is also submodular, and results from Sec- 
tion 3 would remain valid for Z'. The reverse algorithm 
might be better suited to handle problems where an 
optimal number of warehouses tends to be close to the 
total number of potential warehouse sites. 


(17) 


6. Example 


Fig. 3 shows the cost tableau for a 2 X5 X4 multi- 
stage location problem. Fixed charges for warehouses 
1~—5 are: 150, 217, 200,-264 and 140. For an expla- 
nation of how the costs in the tableau are calculated, 
see Fig. 1 and Section 4. The diagram in Fig. 5 is a 
5-dimensional hypercube whose nodes represent all of 
the possible open warehouse combinations. The search 
tree will be a subgraph of this hypercube. The number 
above each node is the value of an optimal solution to 
P when only those warehouses marked in circle of the 
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Demand Centers 


Warehouses 


Dummy 
Factory 


Vig. 3. 


Calculations 

Setting k* = 2 we get from (19) 1332 + 290 
+38 < 1762 < 1332 + 490 + 0. Since /= 

2 = 2 go to Step 4. 


node are open. The number in the parenthesis above a Step 
node is the total cost lower bound for that node 2 
obtained by calculating (16). 

We illustrate the steps of the algorithm applied to 


this example. 4 Let /= 1. Remove {1} from @ (now @ = 9). 

Step Calculations . nee ae on es 3 

1s, Set /= 1, Q =. Z" = Z(Q) = 2107. Applying : 6 PY. 
Rule 1 we get AZ,(/ — {1}) = 2303 — 2201 3 For each i€ {2, 5,3} apply test (2) setting 
= 102; AZ,(/ — {2})=217;AZ3,0 - {I - T= {4} U {2}. All solutions are fathomed 
{3}) = 316; AZ,(/ — {4}) = 162; AZ5(/ — (see numbers in parenthesis in Fig. 5). List(2) 
{5}) = 140.1" =Q=@. For each i€é {3, 2,3, = 9. Go to Step 4. 
4, 5} apply test (16). None of the lower bounds 
exceeds the current upper bound. (See the num- a Let /= 1. Remove (4} from Q (now Q = 9). 
bers in parenthesis in Fig. 5.) Next calculate Let k= 3. Let List(1)= {2, 5}. Let Q= {3}. 
AZ) = 1880 ~ 2107 = A397: AZ,(9) = Let 1 =2.Go to Step 3. 
~94; AZ3(9) = —123; AZ,@ = —136;425@) 3 For each 7€ {2, 5} apply test (16) setting 
= —105. Let List(1) = {2,5,3,4}. Let Q ={1}; T= {3} U {i}. All solutions are fathomed 
Z" = 1880; 1 = 2; go to Step 2. (see Fig. 5). List(2) = @. Go to Step 4. 

2 Z(I) = 2334; Lief; = 1002; the fixed charge 4 Let 2= 1. Remove {3} from Q. (Now Q = @.) 
ordering is f, = 140, f, = 150, f4 = 200, /5 = Let k =5. Let List(1)= {2}. Let Q= {5}. 
264. The DV ordering is DV, = 0, DV, =0, Let ] = 2; Go to Step 3. 
DV, = 38, DV, = 42, DV; = 48. Z" = 1880. : . 
Setting kt =2 os te roe (19) 1332 +290 + 3 For ie {2} apply test (16) setting T= {5} U 
38 < 1880 < 1332 + 490 +0. Since /=1<2, ti}. The | solution is fathomed (see Fig. 5). 
jet / = 2; let Q = {1};g0 to Step 3. List(2) = @. Go to Step 4. 

3 For each i € {2, 5,4, 3} apply test (16) setting : Let != 1. Go to Step 5. 


T= {1} U {i}. Gee Fig. 5.) For each i€ 
{2, 5,4, 3} calculate AZ,({1}). AZ3({1}) = 
118. AZ,({1}) = 32; AZs5({1}) = 22; AZ, 
({1}) = —18. Let Z" = 1762. List(2) = {2}. 
Go to Step 4. 


Stop. Z° = 1762 is optimal for warehouses 
{1,3}. An optimal solution is shown in Fig. 4. 


An optimal solution to the example is given in Fig. 
4 when warehouses 1 and 3 are opened. The greedy 


157 


R.V. Nagethout, G.L. Thompson / Multistage location ~ allocation 


Demand Centers 


Warehouses 


wy 
Y 
od 
u 
° 
w 
is) 
i) 
fy 


Fig. 4. 


2303 


32345 


Fig. 5. 


ba 


solution is optimal and is obtained at the circled ver- 
tex {1,3} in Fig. 5. The upper bound on an optimal 


number of warehouses, k” = 2, was obtained in the 
first application of Step 2 and thus none of level 


three warehouse combinations were examined. In total, 
16 of the vertices in Fig. 5 were generated. Six of the 


vertices on level 2 were fathomed using (16). The 
total number of transportation pivots required to 


solve the sample problem was 85. 


7. Computational results 


The cost operator algorithm of Section 5 was 


coded in FORTRAN IV and the runs were made ona 


Table 1 


Test problems 


Problem Set 


XII 


XTITI 


XIV 


XV 


* Numbers in the brackets correspond to ranges. 


# Problems 


4 


4 


(qx2xr ) 


0x16X50 
0X 16x50 


0X16x50 
0X16x50 
“0x25x50 
0X25x50 
ieee 
0x15x45 
0x15x45 
3X24x50 
3X24x50 
3x24x50 
3X24x50 


3X24x50 


3X24x50 


id. 
il. 


58268 


58268 


58268 
58268 
58268 
58268 
58268 
37440 
37440 
58268 
58268 
58268 
58268 
58268 


58268 


35000 


35000 


35000 


30000 


30600 


30000 


Sosemenesineninnauneaenearenenemermeamemnenes ee sraate ee 
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DEC 20 time sharing system. Many of the problems 
were run at different times during the day and thus 
the execution times may vary up to ten or fifteen per- 
cent depending upon the computing load of the 
machine. The maximum time allowed for solving any 
problem was 500 seconds. 

All of the problems are derived from the Kuehn and 
Hamburger data which was originally presented in [14]. 
Problem sets I through VII are taken from Akinc and 
Khumawala [1],and VIII and IX are taken from Ell- 
wein and Gray [8]. These are single stage problems of 
the type P1. Problem sets X through XV are multi- 
stage problems which are solved here for the first time. 

Since it is necessary to read several articles to deter- 


Ss. f 


—i_ i 
7,590/12,500/ 
5000 17,500/25,000/ 
10000 17,500 
7,500/12,500/ 
15000 17,500/25,000/ 
7,500/12,500/ 
98628 = 17°500/25,000/ 
7,500/12,509/ 
5000 17,500/25,.000/ 
7,500/12,500/ 
15000 17,,500/25,000/ 
7,500/12,500/ 
98628 = 1 7500/25 ,000/ 
1000 15,000 
5000 40,000 
1500 22,500 
7500 60,000 
7,500/12,500/ 
5000 17,500/25,000/ 
7,500/12,500/ 
15000 17,500/25,000/ 
7,500/12,500/ 
38628 = 17 500/25, 000/ 
7,500/12,500/ 
5000 17, 500/25,000/ 
7,500/12,500/ 
15000 17,500/25,000/ 
7,500/12,500/ 
38268 = 17500/25,000/ 
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mine how the problems I through VII were originally 
created, we will describe them here. The Kuehn and 
Hamburger data represents a multistage system with 3 
factories, 24 potential warehouse sites, and 50 demand 
centers located in the continental United States. Prob- 
lems I through IV were formed by considering only 
the first 15 potential warehouse sites (i.e., Atlanta, 
Boston, through New Orleans as in [14]). To obtain 
the shipping costs (c,;) from these sites to the 50 
demand centers multiply each distance by $ 0.025/mile 
unit cost, which is the bulk shipping rate. The sixteenth 
warehouse is actually the factory at Indianapolis. We 
used the distance from Indianapolis to the 50 dernand 
centers multiplied by $ 0.0125/mile unit as its shipping 
cost. The factory at Indianapolis has the same capacity 
as the other warehouses and it has a zero fixed charge. 
The problems V through VII are set up in a similar 
fashion except that all 24 potential warehouse sites are 
used. Problems X through XV are multistage problems 
which are also derived from the Kuehn and Hamburger 
data. All three factory sites were used. The factories 
have capacities of 30 000 or 35 000 units. To derive 


Table 2 
Computational results 


# Vertices 


h to Total Total Total A&k Nauss 

Problem error find # * t # Né&T ‘IBM IBM, 

Set Greedy Optimal Vertices [I | Warehouse Pivots (DEC-20) (570/165/ (azoyrea/ 
I-l 0 25 25 11 13 163 1.61 10.21 7.6 
I-2 0 21 21 11 12 152 1.43 9.15 1.5 
I-3 0 21 21 11 12 156 1.46 9.26 6.8 
I-4 0 21 21 11 12 156 1.45 9.58 10.8 
I 19 58 147 5 8 1712 6.16 19.68 6.6 os 
ILI-1 0 28 32 9 ‘1 158 1.02 23 1.5 
ILI-2 0 30 37 7 9 217 1.19 -43 1.7 
III-3 0 51 102 4 7 858 3.37 38.65 5.3 
III-4 0 41 158 3 5 2582 8.36 34.39 2.9 
Iv-1l 0 28 31 9 11 115 96 47.5 - 
Iv-2 0 27 31 7 9 149 1.03 -44 - 
IV-3 0 38 41 3 5 218 1.2 .30 - 
IV-4 0 30 30 3 4 133 -85 15 - 
y-1 47 254 1303 ll 17 11700 56.7 23 18.4 
V-2 ete 128 6278 9 14 96462 422.8 120* 18.1 
V-3 - - - - - - 500* 120* 9.1 
V-4 - - - - - - 500* 120* 73.4 
VI-1 1G 67 360 10 15 141? 10.35 75 2.2 
VI-2 0 78 908 6 11 6952 31.8 7.75 9.2 
VI-3 19 982 2172 5 8 19420 79.9 120* 18.1 
VI-4 15 1299 16367 2 7 146900 495.2 120% 27.2 
VII-1 0 62 238 10 14 978 5.5 68 - 
VII-2 0 69 408 6 11 2408 9.87 1.65 - 
VII-3 20 104 264 4 8 1922 6.72 1.34 - 
VII-4 06 75 lol 2 4 1088 4.05 -46 Cees 
Vill Q 60 804 5 11 6247 _ 29.6 42.55 1.8 
IX 3.7 92 1304 2 7 9367 41.6 4.17 8.1 


: Computation terminated 
— Not available 
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the shipping costs, we multiplied the distances from 
each factory to each warehouse by $ 0.0125/mile unit 
cost, and from each factory to each demand center by 
$ 0.0375/mile unit cost. To get the shipping costs from 
the warehouses to the demand centers, we multiplied 
the distance by $ 0.0250/mile unit cost. 

Table 1 contains a description of all of the test prob- 
lems. Table 2 contains the computational results for 
the single stage problems and Table 3 for the multi- 
stage problems. The percent error formula used to eva- 
luate the greedy solution was, 


ZO 78 
z* 


where Z° is the greedy value and Z* is the optimal 
value. Also included in Table 2 are the execution times 
for the same problems solved in [1,17]. For test pur- 
poses the problem sets I through XV are very interes- 
ting because they contain problems with a wide range 
of difficulty. 

The results for the greedy solution were amazing. 


60 R.V. Nagelhout, G.L. Thompson / Multistage location — allocation 


jaa 


Table 3 
Computational results 


‘to # Vertices Total Optimal Total CPU 
Problem error to # a # # Time 
Set Greedy Optimal Vertices lr | Warehouses Pivots (Dec-20) 
X-1 0 56 143 11 15 2414 11.88 
X-2 .30 514 858 8 13 12287 61.39 
X-3 24 277 680 6 11 7411 39.15 
X-4 1.6 492 3627 1 7 39084 197.24 
XI-1 0 54 99 9 13 1708 8.51 
XI-2 ) 71 52 5 8 7895 36.09 
XI-3 0 62 184 5 8 2830 13.69 
XI-4 0 90 837 2 6 14523 66.88 
XII-1 .02 66 76 9 12 2874 17.69 
XI1-2 0 69 352 5 8 6885 34.43 
XII-3 ) 75 452 4 8 8431 41.26 
XII-4 0 81 499 3 6 10383 50.1 _ 
XI1I-1 0 56 145 11 14 2412 12.16 
XIIT-2 .30 514 858 8 13 12144 62.97 
XIII-3 124 277 680 6 11 7840 45.44 
XIII- 1.7 492 3627 1 7 41054 246.83 
XIV-1 0 54 100 9 13 2060 12.85 
XIV-2 0 71 537 5 8 8480 44.7 
XIV-3 0 62 184 5 8 3152 19.9 
XIV-4 0 90 838 2 6 14463 75.7 
XV-1 .01 73 84 9 12 3627 23.8 
XV-2 ) 69 363 5 8 7813 44.3 
XV-3 0 75 454 4 8 9187 49.8 
Xv-4 0 81 499 2 6 11234 60.53 
The largest percent error incurred was 3.7 (see Table by three different sets of authors and are difficult to 
2, [X), however in general the percent was less compare due to the differences in computers and in 
than 0.5 percent and the greedy solution was an opti- programming efficiency. As far as machine speeds is 
mal solution in 32 out of 51 problems. concerned, the IBM 370/168 is the fastest, followed 
Another interesting statistic is the number of ver- by the IBM 370/165 and then by the DEC-20. How- 
tices required to find the optimal solution as compared ever, actual speed ratio factors for pairs of the compu- 
to the total number of vertices searched. An optimal ting machines are virtually impossible to find. It seems 
solution is usually located very early in the search fair to state that the performance of the three codes 
process and most of the effort is typically spent veri- shown in the table are not significantly different; we 
fying its optimality. This experience is similar to that may say they represent the state of the art of compu- 
found in solving other kinds of integer programming tational results on these problems. 
problems. Table 3 exhibits computational results on multi- 
The set /* , which represents the warehouses fixed stage problems having 3 factories, 24 warehouses, and 
open in the initialization phase of the cost operator 50 demand centers. As can be noted, the computation 
algorithm, often accounts for as much as 85% of the times are relatively small and exhibit a relatively small 
total number of open warehouses in an optimal solu- variance for integer programming problems. We again 
tion. The uncapacitated problems (sets IV, VII, XII, found the performance of the greedy solution to be 
XV) were very easy to solve as expected. As far as even better than for the single stage problems. This is 
problem difficulty is concerned, those problems for perhaps due to the fact that the factories have a large 
which an optimal number of warehouses is approxi- total capacity, and they can ship to customers directly 
mately one half of the total number of potential if many of the warehouses are closed down. 


warehouse sites are usually the most difficult. This is 
probably due to the fact that the number of possible 


warehouse combinations, ( ,) when / is even or 8. Conclusions 
( ne j.) When J is odd, is the maximum of the bino- 
mial coefficients (/). We have presented a cost operator algorithm for 


The CPU times quoted in Table 2 were obtained solving multistage location-allocation problems which 


R.V. Nagelhout, G.L. Thompson / Multistage location — allocation 161 


does not employ problem relaxations as do the other 
currently best approaches [1,17]. Computational 
results indicate that this method is competitive with 
the others. The greedy solutions obtained by the 
method are usually extremely close to or are optimal. 
Also, because the method computes maay near opti- 
mal solutions as it solves the problem, these near opti- 
mal solutions can be saved and printed out for use by 
a manager if he or she desires. Computational results 
on the solution of multistage location problems are 
presented here, but we have been unable io find other 
published results on such problems for comparison. 
The performance of our method on these problems 

is encouraging. 
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